Surface tension fluctuations and a new spinodal point in glass-forming liquids 
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The dramatic slowdown of glass- forming liquids 
[1] has been variously linked to increasing dy- 
namic [2, 3] and static [4, 5] correlation lengths. 
Yet, empirical evidence is insufficient to decide 
among competing theories [6—11]. The random 
first-order theory (RFOT) [7] links the dynamic 
slowdown to the growth of amorphous static or- 
der, whose range depends on a balance between 
configurational entropy and surface tension. This 
last quantity is expected to vanish when the tem- 
perature surpasses a spinodal point beyond which 
there are no metastable states. Here we mea- 
sure for the first time the surface tension in 
a model glass-former, and find that it vanishes 
at the energy separating minima from saddles, 
demonstrating the existence of a spinodal point 
for amorphous metastable order. Moreover, the 
fluctuations of surface tension become smaller for 
lower temperatures, in quantitative agreement 
with recent theoretical speculation [5] that spa- 
tial correlations in glassy systems relax nonexpo- 
nentially because of the narrowing of the surface 
tension distribution. 

It is only recently that a static correlation length £ 
has been clearly detected in a glassforming liquid [4, 5], 
by measuring a point-to-set correlation function [12, 13]. 
The idea is to consider a liquid region of size R subject 
to amorphous boundary conditions provided by the sur- 
rounding liquid frozen into an equilibrium configuration. 
The external particles act as a pinning field favouring 
internal configurations which best match the frozen ex- 
terior. Clearly, the effect of the border on the innermost 
part of the region is smaller as R grows larger. Less 
trivially, on lowering the temperature the effect of the 
amorphous boundary conditions propagates deeper into 
the region. More precisely, if we measure some correla- 
tion (or overlap) q c {R) between the initial configuration 
of the region and that reached at infinite time under the 
effect of the amorphous boundary conditions, the decay 
of q c (R) is slower the lower T [4, 5]. This demonstrates 
the existence of an increasing static correlation length £. 
Irrespective of its precise definition, q c {R) is the point- 
to-set correlation function. 



According to RFOT, the decay of q c (R) is regulated 
by competition between a surface energy cost, YR e , try- 
ing to keep the region to the same amorphous state as 
the external configuration, and a configurational entropy 
gain, THR , favouring a transition to another of the 
exponentially many states, Af(R) ~ cxp(i? d S), avail- 
able to the region [4]. The surface tension Y is pre- 
dicted to vanish for temperatures higher than a spinodal 
value. Above this point metastable states merge into 
a single ergodic state. The cost/gain terms balance at 
R = £ = (Y/TE)*^: for R < £ the surface cost Y c R e 
keeps the region in the same state as the external envi- 
ronment; for R > £ the entropic gain TY>R d dominates 
and the region is free to rearrange into some other state. 
The RFOT length £ is not only the typical size of the 
rearranging regions; it also represents the largest scale 
over which it is sensible to define a metastable state: a 
state defined over a region much larger than £ is unstable 
against fragmentation into sub-regions of typical scale £. 
For this reason RFOT is also called mosaic theory. 

The existence of many metastable states is central to 
RFOT. But for more than one state to exist, a nonzero 
surface tension is necessarily required. In fact, one may 
argue that the surface tension is as much a fundamental 
ingredient of the theory as the configurational entropy, if 
not more. Hence a direct empirical investigation of the 
surface tension is needed in order to put a more stringent 
test on RFOT. This is what we do here. 

Surface tension is defined as the free energy cost per 
unit area associated to a surface separating two phases 
[14]. The determination of the surface free energy be- 
tween metastable states is very challenging, first, be- 
cause metastable states are amorphous and the interfaces 
are hard to detect, and second, because their lifetime 
is necessarily finite. Excitations are constantly form- 
ing and relaxing: this is the relaxation mechanism of 
RFOT, through which the whole liquid state is slowly 
explored. For our initial study of surface tension we 
choose to avoid such complications and focus on inherent 
structures (ISs), i.e. local minima of the potential energy 
[15]. This approach measures energy rather than free en- 
ergy, so we expect our results to apply quantitatively only 
at quite low temperatures. 



2 





100 





FIG. 1: Map of the excitation. Overlap, i.e. similarity, be- 
tween the configuration right after the exchange of the spheres 
and the hybrid inherent structure after minimization. Only 
a thin slice located at half height of the whole three dimen- 
sional system is shown. Each particle is coloured according to 
how much it moved after the artificial excitation was created 
(dark blue: small displacement, high overlap — ligth blue: 
large displacement, small overlap). Upper panel: two config- 
urations at high temperature T = 1.33 T c . Lower panel: two 
configurations at low temperature T = 0.89 T c . The hybrid 
minimum clearly bears memory (high overlap, dark blue) of 
the parent ISs far from the boundary of the sphere (white 
circle). On the other hand, the memory is lower (low over- 
lap, light blue) along the interface, where particles have been 
moved the most by the minimization procedure. Although 
at lower temperature the interface is somewhat sharper, it 
is in general quite rough. Strong fluctuations of the overlap 
along the interface are evident, clearly indicating that there 
are large surface tension fluctuations. Note that the hybrid 
configuration is nevertheless a typical IS of the system: the 
only reason why we can visualize the interface is that we know 
a priori the shape and position of the excitation and we can 
use the parent configurations as reference to calculate overlap. 
Without this information, it would be impossible to distiguish 
the hybrid inherent structure from any other one. 



We consider a set of ISs obtained minimizing equi- 
librated instantaneous configurations of a model soft- 
sphere glassformer (sec Methods). Our lowest working 
temperature is T = 0.89 T c (T c is the Mode Coupling 
temperature [16]) 

Given a pair a and j3 of ISs, we exchange between them 
all particles located within a sphere of radius R. We 
then minimize the two configurations thus obtained to 
produce two new ISs. Each of them is a hybrid minimum, 
resembling the "parent" minima far from the surface of 
the sphere but rearranged close to it (Fig. 1). 

For each hybrid IS we compute the surface energy E s a g 
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FIG. 2: Surface energy. E s as a function of the excitation 
radius 7? at various temperatures. A single power law fit R e 
would give an exponent 9 > 2, which is hard to justify. On 
the other hand, especially at high T it is clear that there are 
small R corrections to a single power law. Full lines are fits 
using equation (1), with 9 = 2 and u) = 1.5. Dotted lines 
show as a reference the power laws with 9 = 2 and 9 = 3/2. 
Error bars are smaller than symbol size. 



by subtracting the parent ISs' contribution (see Meth- 
ods). Fig. 2 shows the sample-averaged surface energy, 
E s , vs. R for several temperatures. There is a well- 
defined relationship between surface energy and size, reg- 
ulated by the temperature: at fixed R, E s increases by 
decreasing T. The data do not correspond to a single 
power-law scaling, so we propose 



E s = Y no R > - SR U 



(1) 



where is the asymptotic surface tension and u < 8. 
The subleading 5R U correction is quite natural. It is 
present in liquids (with uj = 1) due to curvature effects 
[14], and in disordered systems, where it may arise ei- 
ther from bulk contributions, as in the random field Ising 
model [17], or from interface roughening that lowers the 
surface energy, as in the random bond Ising model [18]. 

We must choose the exponents of eq. (1) with some 
criterion, because the nonlinear fit with four parameters 
is marginally stable, and many sets of parameters give 
good fits. Our data strongly suggest the conservative 
choice 8 = 2, which seems to describe the large R be- 
haviour better than the alternative 8 = 3/2 predicted by 
a wetting argument [19, 20] . The value 8 = 2 is also 
found in spin models with finite range interactions [21]. 
To fix w, we take eq. (1) as valid for the whole popula- 
tion of surface energies (instead of just the average), and 
ascribe all fluctuations of E s a p to the quantity 



YooR 2 - E\ 



a/3 



R u 



(2) 



We then require that the variance of S a p be independent 
of R, which is the typical behaviour of random systems 
[18, 22]. This procedure (see Methods) gives u = 1.5(2). 
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FIG. 3: Surface tension vs. temperature. Upper panel: 
Yoo as a function of the quenching temperature of the inherent 
structures. The vertical dotted line marks the mode coupling 
temperature. The surface tension decreases on increasing T, 
although too smoothly to indicate a sharp spinodal tempera- 
ture. Lower panel: inherent structure energy as a function of 
the quenching temperature. A comparison of the two curves 
clearly shows that they are correlated. 



With and w thus fixed, eq. (1) fits the E S (R) data 
very well (Fig. 2), and we obtain the asymptotic surface 
tension as a function of T (Fig. 3, top). We find 
that Foo decreases for higher T, and becomes quite small 
above T c . Yet, the decrease is rather smooth, so that it 
is hard to define a spinodal temperature. 

On the other hand, the plot of the average IS energy, 
eis, vs. T (Fig. 3, bottom) suggests using the energy as 
a control parameter. Indeed, the vs. eis(T) curve is 
nearly linear (Fig. 4, left), indicating that Y^ vanishes 
quite sharply at a well-defined energy e t h- Interestingly 
enough, eth is very close to the threshold, i.e. the energy 
below which minima start to dominate the energy land- 
scape [23, 24]. The threshold e t h is defined as the point 
where the instability index of saddles vanishes (Fig. 4, 
right). Hence e t h is the true spinodal point of amor- 
phous order, fixing the upper limit of stability of the 
RFOT mechanism. 

We expect the surface tension to fluctuate when con- 
sidering different pairs of ISs. We thus compute the dis- 
tribution of the single-sample surface tension, 
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FIG. 4: The spinodal point. Left: Yoo (squares) vs. IS en- 
ergy. Right: intensive saddle instability index k (circles) vs. 
IS energy. Lines are linear fits to the data. Both the surface 
tension and the instability index seem to vanish at a simi- 
lar energy, the threshold e t h, which is therefore the spinodal 
point. Instability index data are from ref. 24. 
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FIG. 5: Surface tension distribution. Normalized his- 
tograms of the surface tension for T = 0.89 T c at two values 
of R. Note that this distribution does not depend on the expo- 
nent ui and that its average tends to YJx>(T) for R — » oo. The 
distribution narrows on increasing R, and thus on lowering the 
temperature. Solid lines correspond to fits with the Weibull 
distribution P(Y) = -d/dY exp[-(F/?/ c ) e ] (see Ref. 5). In- 
set: variance of Y a p vs. R. For small R this quantity seems to 
decay as 1 / i? 2 , but at larger R it is clearly larger. Note that 
that, Var[Y a/3 ] = Var[Y«,] + Var[<5 Q(3 ]/# 4 ~ 2 "- If Var[<M does 
not depend on R and u> ~ 3/2, then Var[<5 a( 3]/7? 4-2 " ~ 1/R, 
which is compatible with the large R behaviour of the data. 



The distribution P(Y, R, T) for two values of R is shown 
in Fig. 5. The first thing we notice is that the distribu- 
tion is quite broad. This means that the original RFOT, 
which assumed a sharp value of Y, cannot hold strictly. 
If there is a single surface tension, a region smaller than 
£ cannot rearrange, so that a nonfluctuating surface ten- 
sion implies a sharp drop of the point-to-set correlation 
q c (R) at R ~ £. With a fluctuating surface tension, on 



the other hand, any region can decorrelate, as long as 
there are target states with surface tension Y < TY,R d ~ e . 
Our finding is thus consistent with the smooth decay of 
the correlation observed numerically [4, 5]. 

The second important result is that the distribution 
narrows as R grows. To interpret this correctly, note that 
P(Y) depends on T both directly through the quench 
temperature of the ISs and indirectly through R. Since 
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the only relevant lengthscale is £(T), the typical size of 
the rearranging regions, we can define 



Carlo with particle swaps [28]. The pair potential is 



P(Y,T) = P(Y,R = Z(T),T), 



(4) 



where £(T) is taken from simulations [5]. When T de- 
creases, £(T) grows, P(Y,T) narrows, and the decay of 
q c (R,T) sharpens. This is indeed the numerical finding 
[5]: the decay is well described by a compressed expo- 
nential, 



q c (R,T)~exp -(R/Z) 



(5) 



where the anomaly ( grows from 1 at high temperature 
up to ~ 4 at low T (higher £ means sharper decay). The 
link between P(Y, T) and sharpness of decay can be made 
more quantitative. The overlap q c {R) is given by [5] 



q c (R,T)= / P(Y,T;R)dY. 



(6) 



Using the approximation (4) (see Methods) and defining 
y = TY>R d - e we obtain 
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q c (y,T)= / P(Y,T)dY. 



(7) 



In this way we avoid the complexity S, whose esti- 
mate is always very tricky [25]. If eq. (5) holds, we 
must have a similar compressed exponential form, q c ~ 
exp[— (y/yc) ]^ with the two anomalies related by £ — 
Q/{d — 9). We computed q c (y,T) at the lowest temper- 
ature, T = 0.89 T c , where the the difference between en- 
ergy and free energy should be less harmful. A com- 
pressed exponential fit gives ( w 4.2(1), encouragingly 
close to the value £ = 4.0(6) of ref. 5. The agreement 
on the values of £ found with two completely different 
protocols puts the generalized RFOT on a firmer basis. 

The present study provides independent evidence of 
the existence of a surface tension distribution with the 
properties required by a generalized RFOT, and it thus 
supports its validity as a description of deeply super- 
cooled liquid. The fact that P(Y,T) broadens and clus- 
ters around Y ~ for T > T c , (i.e. for E > E th ) is 
significant in two ways. First, it implies a crossover from 
many states to one state through a spinodal mechanism: 
a null surface energy cost means that rearrangements are 
no longer excitations. Secondly, it says that this transi- 
tion is smooth: instead of disappearing abruptly at T c (as 
in mean-field [26]) metastable states slowly "fade out", 
because excitations are becoming less and less costly on 
average and because states are slowly merging with each 
other, as indicated by many pairs of states having near- 
zero surface tension. 



Methods 

a. System. We have simulated a three dimensional 
soft-sphere binary mixture [27] using Metropolis Monte 



Vij(n,rj) = [(<Ti + <Tj)/\Yi 



with a polynomial 



smooth cut-off at large distances (parameters as given 
in Ref. 5). The mode-coupling temperature for this sys- 
tem is T c =0.226 [29]. A system of N = 16384 particles 
was considered in a box of lenght 25.4 (in such units the 
number density is 1). 

b. Minimization. Minimization of instantaneous 
configurations was done with the LBFGS algorithm [30] . 
After exchanging the particles of the sphere, some pairs 
of particles end up at very short distances, giving very 
high energies and gradients that tend to destabilize the 
minimizer. For these configurations, we have minimized 
with a combination of 100 standard Metropolis Monte 
Carlo steps at T = plus LBFGS. The sphere is placed 
so that the density and composition of the resulting con- 
figurations are identical to those of the initial ISs. 

c. Definition of surface energy. We indicate with 
E a fj and Ep a the energies of the hybrid configurations. 
The first label indicates the original parent configuration 
outside the sphere. The surface energy is 
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(8) 
(9) 



d. The exponent oj. Our estimation of oj relies on 
the physical hypothesis that the prefactor S in the cor- 
rection term of Eq. (1) (unlike Y^) has large fluctuations 
even for R — > oo. This behaviour is rather ubiquitous in 
disordered media [17, 18, 22] and finds some confirmation 
here by the violation of the central limit theorem shown 
by the variance of Y a p at large R (see Fig. 5, inset). We 
fix oj at the value for which the variance Var[5(w)] is inde- 
pendent ofi? with the following self-consistent procedure: 
for a running value of the exponent, Q, we fit the aver- 
age of E s via Eq. (1) to obtain Y^lo) and a variance 
Var[5(<D)]. At large R 



Var[£ s ] - Ybx[5{uj)]R 



■2Q 



(10) 



where Var[a3] will in general depend on R. A similar 
formula holds for ui = oj, but with i?-indepdent variance, 
so that 

logVar[5(w)] = logVar[5(w)] + 2(lu - u)logR. (11) 

The procedure is to fit logVar[5] vs. logi? for several 
values of uj to obtain a slope a(u>) which is finally fitted 
to a{Cb) = 2(ui — 6b) to obtain oj. When we do this we do 
not find any trend of oj with the temperature. Thus, to 
get rid of unwanted thermal noise in the determination 
of this exponent we use the lowest T to fix its value. This 
procedure gives oj — 1.5(2), which is very much within 
the range of values provided by nonlinear fit of the data 
surface energy with all four parameters free. 

e. Calculation of the anomaly. The distribution 
P(Y; R) is peak-shaped in Y; as R is increased the peak 
narrows and shifts to the right, approaching Y^. Eq. (6) 
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says that the overlap is the area of P(Y; R) to the right of 
TSi?. Around £, q c (R) is rapidly decaying, so that TT,R 
must be greater than the peak position for all R > £, 
and the integral is just measuring the area of the right 
tail of P(Y;R). To study the shape of q c {R) near the 
inflection point £ (which dominates the value of £), we 
may thus approximate P(Y;R) with P(Y;£) (with the 
consequence of slightly overestimating the sharpness of 
the decay). 

/. The minimal energy with frozen environment. We 
have treated the external and the internal parts of the 
spherical excitations symmetrically. However, one may 
argue that when a region rearranges, it chooses the tar- 
get state with the minimum energy with respect to the 
external environment (S. Franz, private communication). 
This is certainly the case in the numerical experiments of 
ref. 5, where the amorphous boundary of the rearranging 
region is really frozen. We therefore define the minimal 
surface energy at fixed external IS, 

Efa, a =rfL{E a() }-E?-E%, (12) 



where /?o is the (3 that minimizes E a p. From this mini- 
mum surface energy at fixed a, we can obtain the distri- 
bution Po(Y) by varying the external IS a. This Po(Y) 
is the distribution that must be plugged into the general- 
ized RFOT expression for the overlap, (6), to get (. The 
value of the anomaly obtained in this way is £ « 4.4, 
comparable to the value £ w 4.2 found with P{Y). In- 
deed, apart form the poorer statistics that makes Po(Y) 
considerably more noisy than P(Y), the form of the two 
distributions is very similar. 
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